{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# BASIC SEMI-SUPERVISED LEARNING MODELS"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sometimes pressing Enter in TensorFlow is not enough to get good results. Sometimes you have too little data to work with and maximum likelihood estimation gives you some bizarre answer. Or even more common - data and labeling are *expensive* and the relationship between the number of instances in your dataset and quality of the final model is not linear at all.\n",
    "That's when semi-supervised learning (SSL) come to play.\n",
    " \n",
    " > Semi-supervised learning is a class of machine learning tasks and techniques that also make use of unlabeled data for training – typically a small amount of labeled data with a large amount of unlabeled data [[source]](https://en.wikipedia.org/wiki/Semi-supervised_learning).\n",
    " \n",
    "These types of ML algorithms are using both labeled and unlabeled data to get more accurate predictions (or in rare cases to get them at all).  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<center><p><img src=\"https://s3-ap-south-1.amazonaws.com/av-blog-media/wp-content/uploads/2017/09/20182516/dataiku-hadoop-summit-semisupervised-learning-with-hadoop-for-understanding-user-web-behaviours-12-638.jpg\" title=\"Types of learning\"/></p></center>  \n",
    "\n",
    "Img 1 - Types of learning [[source]](https://www.slideshare.net/Dataiku/dataiku-hadoop-summit-semisupervised-learning-with-hadoop-for-understanding-user-web-behaviours)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And since we are dealing with unlabeled data (obviously, unknown instances) we grab some assumptions from unsupervised learning:\n",
    " \n",
    " - Continuity assumption (*Points which are close to each other are more likely to share a label*)\n",
    " - Cluster assumption (*The data tend to form discrete clusters, and points in the same cluster are more likely to share a label*)\n",
    " - Manifold assumption (*The data lie approximately on a manifold of much lower dimension than the input space*)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "They seem reasonable enough though it's not a problem to find counterarguments.\n",
    "Why do we need all these complications? Well, because it's working (in theory)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<center><p><img src=\"https://ai2-s2-public.s3.amazonaws.com/figures/2017-08-08/19bb0dce99466077e9bc5a2ad4941607fc28b40c/4-Figure1-1.png\" title=\"SSL explanation\"/></p></center>\n",
    "    \n",
    "Img 2 - Semi-supervised explanation [[source]](http://people.cs.uchicago.edu/~niyogi/papersps/BNSJMLR.pdf)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The whole idea of SSL is kinda bayesian in its nature, but we won't dive deep.\n",
    "\n",
    "As analogy to SSL [some resources](https://deepai.org/machine-learning-glossary-and-terms/semi-supervised-learning) give human inductive logic: we know just some instances of the phenomenon, but we try to infer general principles and same here: algorithms are using few reference points (our labeled data) to find the general pattern that will suit our unlabeled data best. But it's not always the case. Some algorithms will try to find generalization for the data given (infer a function, that describes our data best), but some algorithms will use thing called [transductive learning](https://en.wikipedia.org/wiki/Transduction_(machine_learning)) ([more here](http://www.cs.cornell.edu/courses/cs4780/2009fa/lecture/13-transduction.pdf)). We will use models that take this approach to solve the classification task. \n",
    "\n",
    "All in all inductive learning is trying to get generalization from the given data and only then predict new data. Transductive learning will just try to predict new data given training data, skipping generalization part.\n",
    "\n",
    "Let's code a bit. I will add comments and explain what I'm doing. I guess most of the time it won't be necessary but still."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Basic imports\n",
    "import os\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import seaborn as sns\n",
    "from sklearn.utils import shuffle"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`sklearn` has two SSL algorithms: Label Propagation and Label Spreading. Let's import them.\n",
    "\n",
    "We will also use [`pomegranate`](https://pomegranate.readthedocs.io/en/latest/index.html). It's a bayesian lib with lots of features but we will take only one-two models from it. It's written in so called skl notation - it uses same syntax and methods names as `sklearn`, so you will understand it pretty quickly. To install it run:\n",
    "\n",
    "through pip\n",
    ">pip install pomegranate\n",
    "\n",
    "through Anaconda\n",
    ">conda install pomegranate"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pomegranate as pg\n",
    "from sklearn.datasets import load_breast_cancer\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "from sklearn.metrics import roc_auc_score\n",
    "from sklearn.semi_supervised import LabelPropagation, LabelSpreading"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import warnings\n",
    "\n",
    "warnings.simplefilter(\"ignore\")  # we don't wanna see that\n",
    "np.random.seed(\n",
    "    1\n",
    ")  # i'm locking seed at the begining since we will use some heavy RNG stuff, be aware"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will use this [dataset](https://scikit-learn.org/stable/datasets/index.html#breast-cancer-dataset). Here is `sklearn` [loader](https://scikit-learn.org/stable/modules/generated/sklearn.datasets.load_breast_cancer.html)\n",
    "\n",
    "Atributes:\n",
    "\n",
    "- ID number \n",
    "- Diagnosis (M = malignant, B = benign) \n",
    "\n",
    "Ten real-valued features are computed for each cell nucleus: \n",
    "\n",
    "- radius (mean of distances from center to points on the perimeter) \n",
    "- texture (standard deviation of gray-scale values) \n",
    "- perimeter \n",
    "- area \n",
    "- smoothness (local variation in radius lengths) \n",
    "- compactness (perimeter^2 / area - 1.0) \n",
    "- concavity (severity of concave portions of the contour) \n",
    "- concave points (number of concave portions of the contour) \n",
    "- symmetry \n",
    "- fractal dimension (\"coastline approximation\" - 1)\n",
    "\n",
    "The mean, standard error and \"worst\" or largest (mean of the three largest values) of these features were computed for each image, resulting in 30 features. For instance, field 3 is Mean Radius, field 13 is Radius SE, field 23 is Worst Radius.\n",
    "\n",
    "All feature values are recoded with four significant digits.\n",
    "\n",
    "Missing attribute values: none\n",
    "\n",
    "Class distribution: 357 benign, 212 malignant"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = load_breast_cancer()\n",
    "df = pd.DataFrame(data.data, columns=data.feature_names)\n",
    "df[\"target\"] = data[\"target\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we will just briefly look through our dataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df.info()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df.describe()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will shuffle our dataset since it has order and we don't want that. \n",
    "Also, I will reduce dimensionality by dropping features to simplify everything.\n",
    "\n",
    "Then we will create our `X` and `Y` and split them into three parts: labeled train data (1), unlabeled train data (2) and test data (3). We will drop a lot of features (se and worst columns; area/perimeter is redundant since it's highly correlated to the radius; concave points feature is redundant too and correlates to concavity and compactness).\n",
    "\n",
    "Little warning: result reproduction is only achievable if you shuffle data from its original order (order that it has in csv file). If you try to shuffle from another state (shuffle when you've already shuffled it, etc.) will give you different results."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df = shuffle(df, random_state=1)\n",
    "X = df.drop(\n",
    "    [\n",
    "        \"target\",\n",
    "        \"radius error\",\n",
    "        \"texture error\",\n",
    "        \"perimeter error\",\n",
    "        \"area error\",\n",
    "        \"smoothness error\",\n",
    "        \"compactness error\",\n",
    "        \"concavity error\",\n",
    "        \"concave points error\",\n",
    "        \"symmetry error\",\n",
    "        \"fractal dimension error\",\n",
    "        \"worst radius\",\n",
    "        \"worst texture\",\n",
    "        \"worst perimeter\",\n",
    "        \"worst area\",\n",
    "        \"worst smoothness\",\n",
    "        \"worst compactness\",\n",
    "        \"worst concavity\",\n",
    "        \"worst concave points\",\n",
    "        \"worst symmetry\",\n",
    "        \"worst fractal dimension\",\n",
    "        \"mean area\",\n",
    "        \"mean perimeter\",\n",
    "        \"mean concave points\",\n",
    "    ],\n",
    "    axis=1,\n",
    ")\n",
    "y = df[\"target\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since we have only seven features now we can do something called [pairplot](https://seaborn.pydata.org/generated/seaborn.pairplot.html): it will plot pointplots between each features and their distribution on diagonal. This will help us find correlations and see how we can normalize our features.\n",
    "This will take a time. More features and objects = more time, so don't try it on large datasets. If you want you can try it on whole data set `df` (it will take more time and every graph will be smaller), you will see highly correlated features that I dropped above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sns.pairplot(X)\n",
    "# sns.pairplot(df)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we will merge labeled and unlabeled data.\n",
    "The common practice (most libs will need it that way) is to label unlabeled data `-1`. Some libs will only accept `NaN` as label for unlabeled data. Always check documentation and user guides."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_1, X_2, X_3 = np.split(X, [int(0.1 * len(X)), int(0.5 * len(X))])\n",
    "y_1, y_2, y_3 = np.split(y, [int(0.1 * len(y)), int(0.5 * len(y))])\n",
    "y_1_2 = np.concatenate((y_1, y_2.apply(lambda x: -1)))\n",
    "X_1_2 = np.concatenate((X_1, X_2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "index = [\"Algorithm\", \"ROC AUC\"]\n",
    "results = pd.DataFrame(columns=index)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "logreg = LogisticRegression(random_state=1, class_weight=\"balanced\")\n",
    "logreg.fit(X_1, y_1)\n",
    "results = results.append(\n",
    "    pd.Series(\n",
    "        [\"Logistic Regression\", roc_auc_score(y_3, logreg.predict_proba(X_3)[:, 1])],\n",
    "        index=index,\n",
    "    ),\n",
    "    ignore_index=True,\n",
    ")\n",
    "results"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Just for reference LR timings.\n",
    "\n",
    "*Reminder: df is 569x31 dense matrix*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%timeit\n",
    "logreg_test = LogisticRegression(random_state=1, class_weight=\"balanced\")\n",
    "logreg_test.fit(df, y)\n",
    "logreg_test.predict_proba(df);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Label Propagation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Time to use our first SSL model - `LabelPropagation` (LP). It's quite simple. Intuitively, it looks and sounds like a clustering algorithm that just \"propagate\" labels to nearest data points. But we need to go a little bit deeper.\n",
    "\n",
    "Y is our label matrix, it has size of $(l+u)\\times C$, where l - amount of labeled data points, u - amount of unlabeled data points, C - number of classes. So in binary classification task we will get 2 columns. It's kind of list that we will get from `predict_proba` method.\n",
    "\n",
    "The algorithm will create a fully connected graph where nodes are **all** the data points. The edges between nodes will be weights $w_{i,j}$: \n",
    "\n",
    "$$w_{i,j} = exp(\\frac{d_{i,j}^2}{\\sigma^2})$$\n",
    "\n",
    "where *d* - distance function (euclidean in this case, but in general could be any distance function you want), $\\sigma$ - hyperparameter, that control (shrink) weights.\n",
    "\n",
    "Then we build probabilistic transition matrix T: \n",
    "\n",
    "$$T_{i,j} = \\frac{w_{i,j}}{\\sum_{k=1}^{l+u}w_{i,j}}$$\n",
    "\n",
    "T is just a matrix with probability of every data point to be in class C. Our labeled data points has 1.0 probability to be in class C (since we know their classes) and unlabeled data will get their classes from neighbours (we calculated weights earlier, they depend on distance between 2 points).\n",
    "\n",
    "The whole algorithm looks like this:\n",
    "1. Propagate Y <- TY (we \"propagate\" labels from labeled data to unlabeled)\n",
    "2. Row-normalize Y (value of element in row / sum of all elements values in row)\n",
    "3. Clamp the labeled data (we fixinig our labeled data, so our algorithm won't change probability of class, or in another words change the label)\n",
    "4. Repeat from step 1 until Y converges (we recalculate or distances and weights, that will give us different transition matrix that will change our belief in assigned labels, repeat until process converge).\n",
    "\n",
    "In case of `sklearn` implementation we have a choice of weightening function: RBF (see formula for $w_{i,j}$) or KNN ($1(x' \\in kNN(x))$). KNN is faster and gives sparse representation. RBF is harder to compute and store dense transition matrix in memory, but has more options for tuning. Also sklearn RBF implementation instead of dividing by $\\sigma^2$ multiply it by `gamma`, so it should be float, not integer.\n",
    "\n",
    "To learn more you can read [this](http://pages.cs.wisc.edu/~jerryzhu/pub/CMU-CALD-02-107.pdf) and [this](https://scikit-learn.org/stable/modules/label_propagation.html#label-propagation). `sklearn` also has some cool examples of `LabelPropagation`: [SVM vs LP](https://scikit-learn.org/stable/auto_examples/semi_supervised/plot_label_propagation_versus_svm_iris.html#sphx-glr-auto-examples-semi-supervised-plot-label-propagation-versus-svm-iris-py), [LP demo](https://scikit-learn.org/stable/auto_examples/semi_supervised/plot_label_propagation_structure.html#sphx-glr-auto-examples-semi-supervised-plot-label-propagation-structure-py), [LP digits recognition](https://scikit-learn.org/stable/auto_examples/semi_supervised/plot_label_propagation_digits.html#sphx-glr-auto-examples-semi-supervised-plot-label-propagation-digits-py) and [LP digits active learning](https://scikit-learn.org/stable/auto_examples/semi_supervised/plot_label_propagation_digits_active_learning.html#sphx-glr-auto-examples-semi-supervised-plot-label-propagation-digits-active-learning-py)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now I will define a little function that will give us a plot of ROC AUC of our algorithm depending on our kernel and parameter list."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def label_prop_test(kernel, params_list, X_train, X_test, y_train, y_test):\n",
    "    plt.figure(figsize=(20, 10))\n",
    "    n, g = 0, 0\n",
    "    roc_scores = []\n",
    "    if kernel == \"rbf\":\n",
    "        for g in params_list:\n",
    "            lp = LabelPropagation(\n",
    "                kernel=kernel, n_neighbors=n, gamma=g, max_iter=100000, tol=0.0001\n",
    "            )\n",
    "            lp.fit(X_train, y_train)\n",
    "            roc_scores.append(roc_auc_score(y_test, lp.predict_proba(X_test)[:, 1]))\n",
    "    if kernel == \"knn\":\n",
    "        for n in params_list:\n",
    "            lp = LabelPropagation(\n",
    "                kernel=kernel, n_neighbors=n, gamma=g, max_iter=100000, tol=0.0001\n",
    "            )\n",
    "            lp.fit(X_train, y_train)\n",
    "            roc_scores.append(roc_auc_score(y_test, lp.predict_proba(X_test)[:, 1]))\n",
    "    plt.figure(figsize=(16, 8))\n",
    "    plt.plot(params_list, roc_scores)\n",
    "    plt.title(\"Label Propagation ROC AUC with \" + kernel + \" kernel\")\n",
    "    plt.show()\n",
    "    print(\"Best metrics value is at {}\".format(params_list[np.argmax(roc_scores)]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gammas = [9e-6, 1e-5, 2e-5, 3e-5, 4e-5, 5e-5, 6e-5, 7e-5, 8e-5, 9e-5]\n",
    "label_prop_test(\"rbf\", gammas, X_1_2, X_3, y_1_2, y_3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ns = np.arange(50, 60)\n",
    "label_prop_test(\"knn\", ns, X_1_2, X_3, y_1_2, y_3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now you can define your model separately with the best (or whatever) hyperparameter value and check its score"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "lp_rbf = LabelPropagation(kernel=\"rbf\", gamma=9e-6, max_iter=100000, tol=0.0001)\n",
    "lp_rbf.fit(X_1_2, y_1_2)\n",
    "results = results.append(\n",
    "    pd.Series(\n",
    "        [\"Label Propagation RBF\", roc_auc_score(y_3, lp_rbf.predict_proba(X_3)[:, 1])],\n",
    "        index=index,\n",
    "    ),\n",
    "    ignore_index=True,\n",
    ")\n",
    "\n",
    "lp_knn = LabelPropagation(kernel=\"knn\", n_neighbors=53, max_iter=100000, tol=0.0001)\n",
    "lp_knn.fit(X_1_2, y_1_2)\n",
    "results = results.append(\n",
    "    pd.Series(\n",
    "        [\"Label Propagation KNN\", roc_auc_score(y_3, lp_knn.predict_proba(X_3)[:, 1])],\n",
    "        index=index,\n",
    "    ),\n",
    "    ignore_index=True,\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "results"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's talk about time. You can just look at the results or run command on your machine. On dense matrixes it's okeish, but on sparse representation it's a catastrophy (maybe knn kernel is ok). But for this we have [Sparse Label Propagation](https://arxiv.org/abs/1612.01414)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%timeit\n",
    "rbf_lp_test = LabelPropagation(kernel=\"rbf\")\n",
    "rbf_lp_test.fit(df, y)\n",
    "rbf_lp_test.predict_proba(df);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%timeit\n",
    "knn_lp_test = LabelPropagation(kernel=\"knn\")\n",
    "knn_lp_test.fit(df, y)\n",
    "knn_lp_test.predict_proba(df);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Label Spreading"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next one is `LabelSpreading` (LS). The algorithm is very similar to spectral clustering algorithm *normalized cuts algorithm* ([look here](https://en.wikipedia.org/wiki/Spectral_clustering), [here](https://towardsdatascience.com/spectral-clustering-for-beginners-d08b7d25b4d8) and [here](https://papers.nips.cc/paper/2092-on-spectral-clustering-analysis-and-an-algorithm.pdf)).\n",
    "LS will create a affinity matrix (same as LP weights calculation step): \n",
    "\n",
    "$$W_{i,j} = exp(\\frac{-||x_i - x_j||^2}{\\sigma^2})$$\n",
    "\n",
    "for every $i\\neq j$ and $W_{i,j} = 0$ for every $i = j$\n",
    "\n",
    "Then we will construct the matrix (Laplacian):\n",
    "\n",
    "$$S = D^{-1/2}WD^{−1/2}$$\n",
    "\n",
    "where D - diagonal matrix with its *(i,i)*-element equal to the sum of the i-th row of W.\n",
    "\n",
    "These two steps is just a spectral clustering on all our data. Next two is more interesting:\n",
    " \n",
    "Iterate $F(t+1) = αSF(t)+(1−α)Y$ until convergence, where α is a parameter in (0, 1), F(t) - classifying function\n",
    "\n",
    "Let $F^*$ denote the limit of the sequence {F(t)}. Label each point $x_i$ as a label $y_i = argmax≤F^*_{i,j}$\n",
    "\n",
    "During each iteration of the third step each point receives the information from its neighbors (first term), and also retains its initial information (second term). The parameter α specifies the relative amount of the information from its neighbors and its initial label information. Finally, the label of each unlabeled point is set to be the class of which it has received most information during the iteration process.\n",
    "\n",
    "You can read article about Label Spreading [here](http://citeseer.ist.psu.edu/viewdoc/download;jsessionid=19D11D059FCEDAFA443FB135B4065A6A?doi=10.1.1.115.3219&rep=rep1&type=pdf), `sklearn` user docs are [here](https://scikit-learn.org/stable/modules/generated/sklearn.semi_supervised.LabelSpreading.html)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def labels_spread_test(kernel, hyperparam, alphas, X_train, X_test, y_train, y_test):\n",
    "    plt.figure(figsize=(20, 10))\n",
    "    n, g = 0, 0\n",
    "    roc_scores = []\n",
    "    if kernel == \"rbf\":\n",
    "        g = hyperparam\n",
    "    if kernel == \"knn\":\n",
    "        n = hyperparam\n",
    "    for alpha in alphas:\n",
    "        ls = LabelSpreading(\n",
    "            kernel=kernel, n_neighbors=n, gamma=g, alpha=alpha, max_iter=1000, tol=0.001\n",
    "        )\n",
    "        ls.fit(X_train, y_train)\n",
    "        roc_scores.append(roc_auc_score(y_test, ls.predict_proba(X_test)[:, 1]))\n",
    "    plt.figure(figsize=(16, 8))\n",
    "    plt.plot(alphas, roc_scores)\n",
    "    plt.title(\"Label Spreading ROC AUC with \" + kernel + \" kernel\")\n",
    "    plt.show()\n",
    "    print(\"Best metrics value is at {}\".format(alphas[np.argmax(roc_scores)]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "alphas = [0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9]\n",
    "labels_spread_test(\"rbf\", 1e-5, alphas, X_1_2, X_3, y_1_2, y_3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "alphas = [0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09]\n",
    "labels_spread_test(\"knn\", 53, alphas, X_1_2, X_3, y_1_2, y_3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ls_rbf = LabelSpreading(kernel=\"rbf\", gamma=9e-6, alpha=0.6, max_iter=1000, tol=0.001)\n",
    "ls_rbf.fit(X_1_2, y_1_2)\n",
    "results = results.append(\n",
    "    pd.Series(\n",
    "        [\"Label Spreading RBF\", roc_auc_score(y_3, ls_rbf.predict_proba(X_3)[:, 1])],\n",
    "        index=index,\n",
    "    ),\n",
    "    ignore_index=True,\n",
    ")\n",
    "ls_knn = LabelSpreading(\n",
    "    kernel=\"knn\", n_neighbors=53, alpha=0.08, max_iter=1000, tol=0.001\n",
    ")\n",
    "ls_knn.fit(X_1_2, y_1_2)\n",
    "results = results.append(\n",
    "    pd.Series(\n",
    "        [\"Label Spreading KNN\", roc_auc_score(y_3, ls_knn.predict_proba(X_3)[:, 1])],\n",
    "        index=index,\n",
    "    ),\n",
    "    ignore_index=True,\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "results"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And timings. They are the same, as LP."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%timeit\n",
    "knn_ls_test = LabelSpreading(kernel=\"rbf\")\n",
    "knn_ls_test.fit(df, y)\n",
    "knn_ls_test.predict_proba(df);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%timeit\n",
    "knn_ls_test = LabelSpreading(kernel=\"knn\")\n",
    "knn_ls_test.fit(df, y)\n",
    "knn_ls_test.predict_proba(df);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Naive Bayes"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let us move to `pomegranate` lib. It supports SSL for some models: Naive Bayes, Bayes Classificator, Hidden Markov Models.\n",
    "We will use only Naive Bayes (NB), it's plain and simple.\n",
    "\n",
    ">Naive Bayes classifiers are a family of simple \"probabilistic classifiers\" based on applying Bayes' theorem with strong (naive) independence assumptions between the features [[source]](https://en.wikipedia.org/wiki/Naive_Bayes_classifier)\n",
    "\n",
    "<center><p><img src=\"https://wikimedia.org/api/rest_v1/media/math/render/svg/52bd0ca5938da89d7f9bf388dc7edcbd546c118e\" title=\"NB explanation\"/></p></center>\n",
    "\n",
    "We will use [Bayes theorem](https://en.wikipedia.org/wiki/Bayes%27_theorem) for our predictions. So, we need to find $p(c|x)$ (probability of the class c given sample x). To do it we need to calculate $p(c)$ - just a class probability in general, $p(x)$ - our *evidence* and $p(x|c)$. The last one is very hard to compute since in general we need to take into account conditional dependencies in our data.\n",
    "\n",
    "But what if we assume that all our features are conditionally independent given class c (very strong and mostly wrong assumption)? Well, this would make our life **a lot** easier.\n",
    "We can now just calculate $p(c|x)$ as a product of $p(x_i|c)$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<center><p><img src=\"https://www.saedsayad.com/images/Bayes_rule.png\" title=\"NB explanation 2\"/></p></center>\n",
    "\n",
    "Img 3 - Naive Bayes explanation [[source]](https://www.saedsayad.com/images/Bayes_rule.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For every sample we will pick the most probable class (this is known as the maximum a posteriori or MAP decision rule). Then we can write down our optimization task"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<center><p><img src=\"https://cdn-images-1.medium.com/max/800/1*0izKZvefo1pCCVXK3fmLhg.png\" title=\"NB optimization\"/></p></center>\n",
    "\n",
    "Img 4 - Naive Bayes classifier [[source]](https://towardsdatascience.com/introduction-to-naive-bayes-classification-4cffabb1ae54)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Another assumption that Naive Bayes make use of is equality of all our features. It's wrong too, but with it we don't need to set weights to our features. Although some libs give you that possibility.\n",
    "\n",
    "There a lot of articles about NB, I like this [one](https://jakevdp.github.io/PythonDataScienceHandbook/05.05-naive-bayes.html) and [wiki's one](https://en.wikipedia.org/wiki/Naive_Bayes_classifier).\n",
    "`pomegranate` bayes classifier user docs are [here](https://pomegranate.readthedocs.io/en/latest/NaiveBayes.html)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we will initiate our NB model directly from our data. `from_samples` is the method that allow us to do it and set distribution that will model our data.\n",
    "Another way is to directly pre-initialize distributions and its parameters and then just predict samples classes.\n",
    "Here I picked `ExponentialDistribution` but you can pick whatever you want. The list of supported distributions is [here](https://pomegranate.readthedocs.io/en/latest/Distributions.html). The ones you want here are: `ExponentialDistribution`, `NormalDistribution`, `PoissonDistribution`. You can check others but it will give you an error since algorithm won't converge."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "nb = pg.NaiveBayes.from_samples(pg.ExponentialDistribution, X_1_2, y_1_2, verbose=True)\n",
    "roc_auc_score(y_3, nb.predict_proba(X_3)[:, 1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Close to random. FeelsBad. Let's try some crazy things. We will construct our independent component distribution with n-dimensions (in this case 7). We will cheat a bit, since we already now how our data is distributed (see `sns.pairplot` step). My take on this is below"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "d = [\n",
    "    pg.ExponentialDistribution,\n",
    "    pg.PoissonDistribution,\n",
    "    pg.NormalDistribution,\n",
    "    pg.ExponentialDistribution,\n",
    "    pg.ExponentialDistribution,\n",
    "    pg.PoissonDistribution,\n",
    "    pg.NormalDistribution,\n",
    "]\n",
    "nb = pg.NaiveBayes.from_samples(d, X_1_2, y_1_2, verbose=True)\n",
    "results = results.append(\n",
    "    pd.Series(\n",
    "        [\"Naive Bayes ICD Prior\", roc_auc_score(y_3, nb.predict_proba(X_3)[:, 1])],\n",
    "        index=index,\n",
    "    ),\n",
    "    ignore_index=True,\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "results"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It's still bad since our data is correlated and not equal at all. At least our `logreg` doesn't think so.\n",
    "But still NB could be useful in some cases since it's simple, **very** fast and has interpretability. You can use it as a baseline or as \"scouting\" algorithm."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.hist(logreg.coef_.reshape(7, 1));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Last but not least - timings (as I said - very fast, even on sparse data)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%timeit\n",
    "nb_test = pg.NaiveBayes.from_samples(pg.ExponentialDistribution, df, y, verbose=False)\n",
    "nb_test.predict_proba(df);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Bonus"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And the last one, but not the least. Pseudo-Labeling. I won't copy-paste it, you can do it yourself. This is very simple approach. We will just separate our labeled and unlabeled data, train model on labelled. Then we will sample unlabeled data and predict these samples and add them to labeled data as new ground truth. That's it.\n",
    "You can literally use everything with it: regression models or classification models. In fact original it was designed for neural nets, but it is very versatile."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<center><p><img src=\"https://datawhatnow.com/wp-content/uploads/2017/08/pseudo-labeling-683x1024.png\" title=\"Pseudo-Labeling\"/></p></center>\n",
    "\n",
    "Img 5 - Pseudo-Labeling explanation [[source]](https://datawhatnow.com/pseudo-labeling-semi-supervised-learning/)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Pseudo-Labeling described [here](https://datawhatnow.com/pseudo-labeling-semi-supervised-learning/) with code and graphs, [copypaste on another blog](https://www.analyticsvidhya.com/blog/2017/09/pseudo-labelling-semi-supervised-learning-technique/).\n",
    "Original article is [here](http://deeplearning.net/wp-content/uploads/2013/03/pseudo_label_final.pdf)\n",
    "You can copy it below and play with it!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
